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ABSTRACT 

The ARGO-YBJ air sho'wer detector monitored the Crab Nebula gamma ray emission from 2007 November to 
2013 February. The integrated signal, consisting of ^3.3 x 10^ events, reached the statistical significance of 
21.1 standard deviations. The obtained energy spectrum in the energy range 0.3-20 TeV can be described by a 
power law function dN/dE = Iq (E / 2 TeV)““, with a flux normalization Iq = (5.2 ± 0.2) x 10“'^ photons cm“^ 
s“' TeV“' and a = 2.63 ± 0.05, corresponding to an integrated flux above 1 TeV of 1.97 x 10“" photons cm“^ 
s“'. The systematic error is estimated to be less that 30% for the flux normalization and 0.06 for the spectral 
index. Assuming a power law spectrum with an exponential cutoff dN/dE = Iq (E / 2 TeV)““ exp (-E / Ect), the 
lower limit of the cutoff energy E„« is 12 TeV, at 90% confidence level. Our extended dataset allows the study 
of the TeV emission over long timescales. Over five years, the light curve of the Crab Nebula in 200-day bins is 
compatible with a steady emission with a probability of 7.3 x 10“^. A correlated analysis with Fermi-LAT data 
over ^4.5 years using the light curves of the two experiments gives a Pearson correlation coefficient r = 0.56 ± 

0.22. Concerning flux variations on timescales of days, a “blind” search for flares with a duration of 1-15 days 
gives no excess with a significance higher than four standard deviations. The average rate measured by ARGO- 
YBJ during the three most powerful flares detected by Fermi-LAT is 205 ±91 photons day“', consistent with 
the average value of 137 ± 10 day“*. 

Keywords: gamma rays: stars - pulsars: individual (Crab) 


1. INTRODUCTION 

The Crab Nebula is the remnant of a supernova exploded in 
1054 A.D. at a distance of ^2 kpc. It contains a 33 ms pulsar 
that powers a wind of relativistic particles. The interactions 
of these particles with the remnant gas, photons and magnetic 
held produce a non-thermal radiation extending from radio 


waves to TeV gamma rays (iBiihler & Blanfor320l4 and ref¬ 
erences therein). Most of the emission is generally attributed 
to synchrotron radiation of relativistic electrons and positrons. 
The spectral energy distribution (SED) peaks between opti¬ 
cal and X-ray frequencies. A second component arises above 
^400 MeV, interpreted as Inverse Compton (IC) of the same 
electrons scattering off synchrotron photons and CMB pho- 
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tons. 

The Crab Nebula is one of the most luminous sources of 
very high energy (VHE) gamma rays in the sky, and the 
first source to be detected at TeV energies (IWeekes et al.l 
119891) . Thanks to its high flux and apparent stability it 
is considered a reference source in gamma ray astronomy. 
Detected by many experiment s , both Cherenkov t elescopes 
(lAharonian et al.l 2004 l2 006t lAlbe rt et al.| l2008h and air 


shower arrays dAmenomori et alj|2009tlAbdo et al.H201^ . the 

Crab Nebula was often used to check the detector perfor¬ 
mance, including sensitivity, pointing accuracy and angular 
resolution. 

In 2010 September the AGILE satellite unexpectedly de¬ 
tected a strong flare from the direction of the Crab Nebula at 
energies above 100 MeV. It lasted two days, with a maximum 
flux t hree times higher than the a verage value (iTayani et al.l 
1201 Ih . later confirmed by Eermi (I Abdo et ^ l20 1 ll) ~ Erom 
then on, Eermi and AGILE reported some more flares, char¬ 
acterized by a rapid increase and decay of the flux, typically 
lasting a few days. The most impressive occurred in 2011 
April, when the observed flux was ^10 times higher than 
usual (iBuehler et al.ll201^ . The measured SED shows a new 
spectral component emerging during flares, peaking at high 
energies (up to hundreds MeVs in the 2011 April flare), at¬ 
tributed to a synchrotron emission of a population of elec¬ 
trons accelerated up to energies of 10'^ eV. The Eermi-LAT 
data also show that these sharp emission peaks are superim¬ 
posed to long-lasting smoother modulat ions with timescales 
of weeks or months (IStriani et al.ll2013T) . The observed flux 
variations are attributed to the nebula, s ince the pulsar emis - 
sion was found to be stable within 20% (IBuehler et al.ll201^ . 
However the origin of this activity is still unclear. In this 
scenario, observations at higher energies could provide pre¬ 
cious information to understand the mechanisms responsible 
for this behaviour. 

A preliminary analysis of the data recorded by the air 
shower detector ARGO-YBJ during the flares, showed an 
increase of the Crab flux at TeV energies wit h a moder¬ 
ate statistical signifi cance, in 3 out of 4 flares (lAielli et alJ 
|20 1 Ot f^rnettoll20 1 3l) . These results have not been confirmed 
by Imaging Atmospheric Cherenkov Telescopes (lACTs) be¬ 
cause the Moon light hampered the observations. However, 
sporadic and short measurements carried out during the first 
part of the 2010 September flare by MAGIC and VERITAS 
show no evidence for a flux variability (lMariottill2010t lOngl 
I 2 OIOI) . The observations by VERITAS and HESS during a 
flare in 2013 March (when ARGO-YBJ was already switched 
off), report a counting rate consistent w ith the steady flux 
(I Aim et alJl2013HAbr^owski et al .1120131) . 

In this paper we present a detailed analysis of the ARGO- 
YBJ data, carried out with a better reconstruction of the 
shower arrival direction, obtained applying quality cuts on the 
events. The study concerns not only the flaring episodes, but 
the whole Crab Nebula data set, consisting of more than five 
years of observation. The ARGO-YBJ layout and operation 
mode are presented in Section 2, with a particular attention 
to the performance in gamma ray astronomy. In Section 3 the 
analysis technique to extract the gamma ray signal is outlined, 
followed by the results obtained with Crab Nebula data. Sec¬ 
tion 4 reports the energy spectrum evaluation and discusses 
systematic errors. In Section 5 the analysis of the time behav¬ 
ior of the Crab Nebula signal is presented, with a search for 
possible flares and rate variations on different timescales. A 


time correlation analysis with the Eermi-LAT data at energy 
E > 100 MeV during ~4.5 years is also reported. Einally, 
Section 6 contains a summary of the results and concluding 
remarks. 

2. THE ARGO-YBJ EXPERIMENT 

The ARGO-YBJ is a “full coverage” air shower detector lo¬ 
cated at the Yangbajing Cosmic Ray Laboratory (Tibet, RR. 
China, longitude 90.5° East, latitude 30.1° North) at an alti¬ 
tude of 4300 m above sea level, devoted to gamma ray astron¬ 
omy at energies above ~300 GeV and cosmic ray studies at 
energies above ^ 1 TeV. 

During its lifetime, from 2007 November to 2013 Eebru- 
ary, ARGO-YBJ monitored the gamma ray sky with an in- 
tegrated sensitivity ranging from 0.24 to ^1 Crab Units 
(iBartoli et al.l 120131) and studied in detail the emission of 
the most luminous gamma ray sources at energies above 
300 GeV, namely the Crab Nebula MGRO J1908-H0 6 
(IBartoli et al.llM2l) . ElESS J1841-055 (IBartoli et al.ll20131j) . 
the Cygnus Re gion (I Bartoli e t al.l 12012bl 2014 |k and the 
blazars Mrk401 ( IBartoli et al.ll201 ll : lAielli et al.ll201()bl) and 
Mrk501 (IBartoli et al.ll2012d) . 

The detector consists of a ^74 x 78 m^ carpet made of a 
single layer of Resistive Plate Chambers (RPCs) with ~92% 
of active area, surrounded by a partially instrumented (^20%) 
area up to ^100 x 110 m^. The apparatus has a modular 
structure, with the basic data acquisition element being a clus¬ 
ter (5.7 X 7.6 m2), made of 12 RPCs (2.85 x 1.23 m^). Each 
RPC is read by 80 strips of 6.75 x 61.8 cm^ (the spatial pix¬ 
els), logically organized in 10 independent pads of 55.6 x 
61.8 cm2 yvJjjqJj indivi dually acquired an d represent the 
time pixels of the detector (lAielli et alJl2006^ . To extend the 
dynamical range up to PeV energies, each RPC is equipped 
with two large pads (139 x 123 cm2) collect the total charge 
devel oped by the particles hitting the detector (lAielli et alJ 
l2012h . The full experiment is made of 153 clusters (18360 
pads), for a total active surface of ^6600 m2. 

ARGO-YBJ operated in two independent acquisition 
mode s: the shower mode and the scaler mode (lAielli et all 
I 2 OO 8 I) . In this analysis we refer to the data recorded from the 
digital read-out in shower mode. In this mode, an electronic 
logic was implemented to build an inclusive trigger, based on 
a time correlation between the pad signals, depending on their 
relative distance. In this way, all showers with a number of 
fired pads ^pad > Nf;.,g in the central carpet in a time win¬ 
dow of 420 ns generated the trigger. This trigger worked with 
high efficiency down to N>r,v = 20, keeping t he rate of random 
coincidences negligible (lAloisio et al.ll2004l) . 

The time of each fired pad in a window of 2 jtsec around 
the trigger time and its location were recorded. To calibrate 
in time the 18360 pads, a software procedure has been de¬ 
velop ed, based on the Characteristic Plane method (iHe et all 
l2007h that using the secondary particles of large vertical 
showers as calibration beams, iteratively reduces the differ¬ 
ences betwee n the time measur ements and the time fit of the 
shower front (lAielli et al.ll2009h . 

The full detector was in stable data taking with the trigger 
condition N,„g = 20 and an average duty cycle ^86%. The 
trigger rate was ^3.5 kHz with a dead time of 4%. 

The detector performance and capabilities in gamma ray 
astronomy have been studied and improved through Monte 
Carlo simulations describing the shower d evelopment in the 
atmosphere by using the CORSIKA code (iHeck et al.lll998^ 

































































Crab Nebula by ARGO-YBJ 


3 



Figure 1. Sensitivity in the detector field of view. Solid line: zenith angle 
dependence of the sensitivity to a Crab-like source. The sensitivity is nor¬ 
malized to a zenith angle 0 = 0°. Dotted line: zenith angle dependence of 
the gamma ray event rate. Dashed line: zenith angle dependence of the back¬ 
ground rate; both rates are normalized to 0 = 0°. 


and the detector respons e with a code based on the GEANT 
package dGEANTlIlQQSh . 

2 .1. Field of View 

One of the distinctive features of air shower arrays is the 
large held of view (FOV), in principle including the entire 
overhead sky. Gamma ray sources cross the FOV with differ¬ 
ent paths according to their declinations. The sensitivity is not 
uniform in the held of view. Given a photon Hux, the atmo¬ 
spheric absorption reduces the rate of showers for increasing 
zenith angles. The cosmic ray background also decreases but 
more slowly, and the combination of the two rates determines 
the trend of the sensitivity as a function of the zenith angle. 
FiglT] shows the event rate in ARGO-YBJ expected from a 
Crab-like source as a function of the zenith angle 0, normal¬ 
ized to the rate at 6* = 0°, compared to the background rate. In 
the same hgure the dependence of the detector sensitivity on 
0 is also reported. According to simulations, the sensitivity at 
6 = 30° (45°) is reduced by a factor (^10) with respect to 
the sensitivity at 0 = 0°. 

The capability to detect a given source depends on its path 
in the held of view (determined by the source declination), 
and in particular on the amount of time that the source lies 
at different zenith angles. The maximum signihcance is for 
a declination 5max = A, where A = 30.1° is the latitude of the 
detector. Given a Crab-like source, the sensitivity decreases 
by less than 10% for declinations \5-5max\ <10°, while it is re¬ 
duced by a factor ^2 for declinations \5-5max\ ^30°. The dec¬ 
lination dependence is slightly stronger (weaker) for sources 
with softer (harder ) spectra with respect to the Crab Nebula 
(iBartoli et al.ll201^ . 

At the ARGO-YBJ site, the Crab Nebula (declination 6 = 
22.01°) culminates at a zenith angle 0 = 8.1° and lies at zenith 
angles 9 < 45° for 6.6 hours per sidereal day. In general, fol¬ 
lowing a source for a longer time per day increases the signal 
signihcance, because of the increasing statistics, but since the 



Figure 2. Angular resolution for different Npad intervals, according to sim¬ 
ulations. The curves represent the fraction of events beyond the angular dis¬ 
tance d from the source, as a function of d. 

signal to background ratio decreases at large zenith angles, 
there is a maximum zenith angle beyond which the signih¬ 
cance begins to reduce. According to simulations, the maxi¬ 
mum zenith angle for the Crab Nebula is ^45°. 

2 .2. Angular Resolution 

The sensitivity needed to observe a gamma ray source is re¬ 
lated to the angular resolution, which determines the amount 
of cosmic ray background. We evaluate the shower arrival 
direction by htting the shower front with a conical shape cen¬ 
tered on the shower core position, to take into account the 
time delay of secondary particles with respect to a hat front, 
a delay that increases with the distance from the core. We set 
this delay to 0.1 ns m“' (lAielli et al.ll2009[) . 

The high granularity of the detector allows the study of the 
shower prohle in great detail and the accurate determination 
of the core position by htting the lateral density distribution 
with a Nishimura-Kamata-Greisen-like function. According 
to simulations, the core position error depends on the number 
of hit pads N^a^ and on the core distance from the detector 
center. For gamma ray induced showers with a core distance 
less than 50 m, the average core position error is less than 8 
(2)mforNp«rf> 100(1000). 

The point spread function (PSF) also depends on N^a^, and 
for a given 'FIpad value, it worsens as the shower core distance 
from the detector center increases. The angular resolution for 
showers induced by cosmic rays has been checked by study¬ 
ing the Moon shadow, observed by ARGO-YBJ with a statis¬ 
tical signihcance of ^9 standard deviations per month. The 
shape of the shadow cast by the Moon on the cosmic ray Hux 
provides a measurement of the detector PSF. This measure¬ 
ment has been found to be in excellent agreement with expec- 
tations, conhrming th e reliability of the simulation procedure 
(IBartoli et al.ll201 iH) . 

The PSF for gamma ray showers is narrower than the cos¬ 
mic ray one by ^30-40%, due to the better dehned time pro¬ 
hle of the showers. To improve the angular resolution for 
gamma ray astronomy studies, quality cuts have been imple- 
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Log^p Energy (TeV) 


Figure 3. Normalized distribution of the primary gamma ray energy for different intervals, for a Crab-like source. 


Table 1 

Characteristics of Crab Nebula simulated events. 


^ pad 

Dcut’^ 

(m) 

Core position 
error^ (m) 

R39° 

(deg) 

Median energy 
(TeV) 

20-39 

no limits 

37 

1.88 

0.34 

40-59 

no limits 

28 

1.50 

0.53 

60-99 

90 

12 

1.04 

0.79 

100-199 

70 

6.8 

0.70 

1.3 

200-299 

60 

4.2 

0.50 

2.1 

300-499 

60 

3.3 

0.41 

3.1 

500-999 

40 

2.3 

0.32 

4.8 

1000-1999 

30 

1.6 

0.24 

8.1 

> 2000 

30 

1.0 

0.19 

17.7 


Note. — 

“ Maximum distance of the shower core from the detector 
center heyond which the events are rejected. 

Distance between the true and the reconstructed cores con¬ 
taining 68% of the events. 

Angular resolution, defined as the 39% containment radius. 

merited, by rejecting the events with a core distance larger 
than a given value Dcm (depending on Npa^) and with an av¬ 
erage time spread of the part icles with respect t o the fitted 
shower front exceeding 9 ns (iBartoli et al.l 120131) . The val¬ 
ues of Dcut are given in Table [T] The fraction of gamma rays 
passing the selection cuts depends on Npad averaging ^80%, 
whereas the fraction of surviving background events is ^76% 
forNpaj/ <100 and ^50% forN^a^ > 100. The selection also 
acts as a mild gamma/hadron discrimination for events with 

> 100 (the sensitivity increases by a factor ^1.1). 

The arrival directions of the selected showers are also cor¬ 
rected for the systematic error due to the partial sampling of 
the sh ower front when the core is close to the edge of the de¬ 
tector (lEckmannll991l) . This systematic error is related to the 


angle between the vector “shower core-detector center” and 
the shower arrival direction. For events with Npad > 100, for 
which the core position is determined with more accuracy, the 
error can be considerably reduced. 

These selections and corrections shrink the PSF by a factor 
ranging from ~ 1.1 for events with = 20-39, up to ~2, 
for Nparf > 1000. The PSFs obtained by simulating the Crab 
Nebula along its daily path up to 0 = 45° are shown in Figl2] 
for different intervals of 

To describe the PSFs analitically, that for small values of 
Npai/ cannot be simply fitted by a two-dimensional Gaussian 
function, the simulated distributions have been fitted with a 
linear combination of two Gaussians. In general, when the 
PSF is described by a single Gaussian (F(r) = l/(27rcr^) exp 
(-r^/cr^), where r is the angular distance from the source po¬ 
sition), the value of the root mean square a is commonly de¬ 
fined as the “angular resolution”. In this case the fraction of 
events within 1 cr is 39%. For our PSFs, the value of the 39% 
containment radius ranges from 0.19° for > 2000 to 
1.9° for 'Hpad = 20-39. Table [T] reports the values of R^g for 
different intervals, together with the core position error, 
after quality cuts, as obtained by simulating the source during 
the daily path in the ARGO-YBJ field of view. 

2.3. Energy Measurement 

The number of hit pads is the observable related to 
the primary energy that is used to infer the source spectrum. 
In general, the number of particles at ground level is not a 
very accurate estimator of the primary energy of the single 
event, due to the large fluctuations in the shower development 
in the atmosphere. Moreover, for a given shower, the number 
of particles detected in a finite area detector like ARGO-YBJ 
depends on the position of the shower core with respect to 
the detector center; for small showers this is especially poorly 















Crab Nebula by ARGO-YBJ 


5 


20 

15 

10 

5 

0 

80 82 84 86 88 

Rigth Ascension (deg) 

Figure 4. Significance map of the Crab Nebula region, after smoothing and 
background subtraction. 
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Figure 5. Comparison of experimental and simulated PSFs: event rate per 
solid angle as a function of the angular distance from the source position, for 
events with > 20, 100 and 500. The solid lines are the results of the 
simulation. 


determined. 

The relation between and the primary gamma ray en¬ 
ergy of showers surviving the selection cuts is illustrated in 
Figl3 where the corresponding primary energy distributions 


for different intervals are reported, as obtained by sim¬ 
ulating a Crab-like source with a power law spectrum with 
index -2.63. The distributions are broad, with extended over¬ 
lapping regions, spanning over more than one order of magni¬ 
tude for small Npad values. The median energies for different 
Nparf intervals are given in Table[T] They range from 340 GeV 
for events with = 20-39, to ^18 TeV for > 2000. 

Since the variable Np„d does not allow the accurate mea¬ 
surement of the primary energy of the single event, the energy 
spectrum is evaluated by studying the global distribution of 
Nparf. The observed distribution is compared to a set of sim¬ 
ulated ones obtained with different test spectra, to determine 
the spectrum that better reproduces the data. 

3. THE CRAB NEBULA SIGNAL 

The data set used for this analysis contains all the events 
recorded from 2007 November to 2013 February, with Npaj > 
20. The total on-source time is 1.12 x 10"^ hours. 

For each source transit, the events are used to fill a set of 
nine 12°x 12° sky maps centered on the Crab Nebula posi¬ 
tion, with a bin size of 0.1°x0.1° in right ascension and dec¬ 
lination (“event maps”). Each map corresponds to a defined 
Np^rf interval: 20-39, 40-59, 60-99, 100-199, 200-299, 300- 
499, 500-999,1000-1999 and Np„j > 2000. 

To extract the excess of gamma rays, the cosmic ray back¬ 
ground has to be estimated and subtracted . Using the time 
swapping method (lAlexandreas et al.]ll993h . the shower data 
recorded in a time interval At = 2-3 hours are used to evaluate 
the “background maps”, i.e. the expected number of cosmic 
ray events in any location of the map for the given time in¬ 
terval. This method assumes that during the interval At the 
shape of the distribution of the arrival directions of cosmic 
rays in local coordinates does not change, while the overall 
rate could change due to atmospheric and detector effects. 
The value of the time interval At is less than a few hours to 
minimize the systematic effects due to the environmental pa¬ 
rameters variations that could change the distribution of the 
arrival directions. 

The time swapping method is a sort of “simulation” based 
on real data: for each detected event, nj "fake" events (with nf 
= 10) are generated by replacing the original arrival time with 
new ones, randomly selected from an event buffer that spans 
the time At of data taking. By changing the time, the fake 
events maintain the same declination of the original event, 
but have a different right ascension. With these events a new 
sky map (background map) is built, with a statistics ny times 
larger than the event map in order to reduce the fluctuations. 
To avoid the inclusion of the source events in the background 
evaluation, the showers inside a circular region around the 
source (with a radius related to the PSF and depending on 
Nparf) are excluded from the time swapping procedure. A cor¬ 
rection on the number of swaps is appli ed to take into accoun t 
the rejected events in the source region (iFlevsher et al.ll2004l) . 

Event and background maps are then smoothed according 
to the PSF corresponding to each N^a;/ interval. Finally, the 
smoothed background maps are subtracted to the smoothed 
event maps, obtaining the “excess maps”, where for every bin 
the statistical significance of the excess is calculated as: 

Ne-Nb 
\J SN^ + 6Ng 

with Ne = E; n; w,- and Nb = E,- bi w,- / «/. In these expres¬ 
sions n, and bj are the number of events of the i'* bin of the 
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Figure 6. Crab Nebula differential energy spectrum multiplied by mea- 
sured by ARGO-YBJ and other experiments. The thick solid line (red in the 
online version) represents the best fit of the ARGO-YBJ data. The dotted 
lines delimit the 1 sigma error band of the Milagro spectrum. 


event map and background map, respectively, w, is a normal¬ 
ized weight, proportional to the value of the PSF at the angular 
distance of the i'^ bin, and n f is the number of swappings. The 
sum is over all the bins inside a radius R^ax, chosen to contain 
the signal events and depending on the PSF. Since the number 
of events per bin is large, the fluctuations follow the Gaussian 
statistics, hence the errors on Ne and Nb are: 5Ne = y/’Sitiiwf 
and SNb = /n^j. 

The number of gamma ray events from the source is: 

Ne-Nb 

=_ - _ £ _ 

Itt Jq'"“ wir^-rdr 

where w{r) is the weight used in the smoothing procedure 
calculated at the angular distance r from the source position. 

When adding all data, an excess consistent with the Crab 
Nebula position is observed in each of the 9 maps, with a total 
statistical significance of 21.1 standard deviations. 

The number of excess events are ^3.3 x 10^, correspond¬ 
ing to 189 ± 16 day“*, where a “day” means a source tran¬ 
sit. Table |2] gives the signal significance for each map and 
the corresponding event rates measured from the source. For 
comparison the background rates measured inside an angular 
window of 1 ° radius around the source are given. Fig.|4]shows 
the total significance map. 

Finally, the gamma ray signal can be used to check the de¬ 
tector angular resolution since the Crab Nebula angular size 
is small compared to the detector PSF. Fig. |5] shows the dis¬ 
tribution of the arrival directions of the excess showers with 
respect to the source position, for Npad > 20, 100 and 500, 
compared to simulations. The agreement is excellent. 

4. ENERGY SPECTRUM 

The energy spectrum is evaluated by comparing the num¬ 
ber of events detected from the Crab Nebula in the previously 
defined N^a;/ intervals to the expected number given by a sim¬ 


ulation assuming a set of test spectra. We consider the power 
law spectrum: 

— (/o,a)-/o|^2TeV ) ’ 

where the flux normalization lo and slope a are the parameters 
to be estimated with the fitting procedure. 

The fit is made by minimizing the value of x^, evaluated for 
any couple of parameters as: 


X^(k,a) 




where Ni^ and are the number of events detected and 
expected, respectively, in the /* N^a^ interval. 

The obtained best-fit parameters are Iq = (5.2 ± 0.2) x 10“'^ 
photons cm“^ s“' TeV“\ and a = -2.63 ± 0.05, with = 5.8 
for 7 degrees of freedom, corresponding to a p-value p = 0.56. 
The integral flux above 1 TeV is 1.97 x 10“" photons cm“^ 
s“^ The flux at 1 TeV obtained in this work is 7 % higher than 
that r eported in a previous ARGO-YBJ paper (iBartoli et alJ 
l2013h . The difference is due to the correction of the event 
rates applied in this work, to reduce environmental and detec¬ 
tor effects on the trigger rate, as described in Section 5.3. 

Fig. |6] shows the obtained spectrum compared with the re¬ 
sults of other experiments. The energy of each point is the 
gamma ray median energy for the corresponding Npad inter¬ 
val. The values of energies and differential fluxes are given in 
Table|2l The spectrum is consistent with a constant slope from 
^300 GeV to ^20 TeV and agrees rather well with the mea¬ 
surement by HEGRA and MAGIC, whereas the HESS and 
Milagro fluxes are about 20% higher in the ~1-10 TeV en¬ 
ergy range. 

The data are less clear conc erning a possible e nergy cutoff 
at higher energies. MA GIC (lAIbert et alJl2008h and HESS 
(lAharonian et al.l l2006l) show a steepening below 20 TeV, 
while the HEGRA spectrum is har der and continues wit h 
a slight softening up to ^75 TeV (lAharonian et al.l [200411 . 
A possible cutoff is also observed by Milagro at ^30 TeV 
(lAbdo et al.ll2012^ . The limited statistics of our data at high 
energy does not allow to draw any conclusion about the spec¬ 
tral properties above 20 TeV. Selecting events with Npa^/ > 
3000 (whose median energy is 26 TeV assuming a power law 
spectrum with index a = -2.63) the statistical significance of 
the signal is 0.75. 

When fitting the data with a power law spectrum with an 
exponential cutoff: 


yiv) “p (-£/£») 

the obtained p-value is always smaller than without a cutoff, 
for any value of Ecut- Eor Ecut = 14.3 TeV (the best-fit value 
obtained by HESS) the p-value is 0.13. We found that the 
p-value is larger than 10% for any value of Ecut >12 TeV, 
indicating that the presence of a cutoff above ^10 TeV cannot 
be excluded even if our data seems more consistent with a 
pure power law. 


dN 

'dE 


(/q, (XjEciit) — 7o 


4.1. Estimation of Systematic Errors 

The previous results can be affected by systematic errors 
of different origin. In the following we discuss the possible 
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Table 2 

Summary of Crab Nebula data. 


Nparf Photon rate Background rate" Significance ^med Differential flux 

(events day“*) (events day^') (s.d.) (TeV) (ph cm“^ s“'TeV^*) 


20-39 

56.7 ± 12.2 

1.3x10'* 

4.6 

0.34 

(6.23 ± 1.34)xl0-“’ 

40-59 

75.4 ± 9.2 

1.1x10“* 

8.2 

0.53 

(1.80±0.21)xl0““’ 

60-99 

34.7 ± 3.9 

4.2x10^ 

9.0 

0.79 

(5.92 ± 0.66)X10-“ 

100-199 

15.4 ± 1.7 

1.9x103 

8.9 

1.30 

(1.37 ±0.15) X10-“ 

200-299 

5.23 ±0.61 

4.9x10^ 

8.5 

2.1 

(5.30 ± 0.63)xl0-‘3 

300-499 

3.51 ±0.44 

3.8x103 

8.0 

3.1 

(1.75±0.22)x10-‘3 

500-999 

2.07 ± 0.27 

2.4x103 

7.6 

4.8 

(5.62 ± 0.74)xl0-‘3 

1000-1999 

0.50 ±0.13 

87.4 

3.8 

8.1 

(1.00 ± 0.26) xl0-‘3 

> 2000 

0.23 ± 0.07 

34.2 

3.5 

17.7 

(1.87±0.54)xl0-‘'' 


Note. — 

“ Average background rate within an angular distance of 1° from the source. 


sources of systematics, evaluating their effects both on the 
flux normalization and the spectral slope. 

1) Energy scale. In our measurement the number of hit pads 
Nparf is used as an estimator of the primary energy. The rela¬ 
tion between the primary energy and Npad is given by Monte 
Carlo simulations. Possible uncertainties and simplihcations 
in the simulation procedure (both in the shower development 
and the detector response) could produce an incorrect Npad 
value and consequently an error in the energy scale. 

The energy scale reliability has been checked using the 
Moon shadow. Due to the geomagnetic held, cosmic rays 
are dehected according to their energy and the shadow that 
the Moon casts on the cosmic ray hux is shifted with respect 
to the Moon position by an amount depending on the energy. 
The westward shift of the shadow has been measured for dif¬ 
ferent Npad intervals and compared to simulations. From the 
analysis of the Moon data, we found that the total absolute 
energy scale error is less than 13% in the proton energy range 
^1-30 TeV (iBartoli et al.ll201 fbh . This estimate includes the 
uncertainties of the cosmic ray elemental composition and the 
hadronic interaction model. 

From this result, given a gamma ray spectrum with index 
a = -2.63, the corresponding systematic error in the flux nor¬ 
malization would be less than 22%. 

2) Pointing error. Fitting the angular distribution of gamma 
rays around the Crab Nebula position we found that the point¬ 
ing error is less than 0.1°. A pointing error affects the mea¬ 
sured gamma ray flux, since the number of photons is ob¬ 
tained by a smoothing procedure weighting the events with 
a PSF centered at the source nominal position. An incor¬ 
rect position would produce a loss of signal. Since the PSF 
is narrower for events with large N^a^/, the loss is larger at 
high multiplicities, and generates a steepening of the spec¬ 
trum. According to our simulation, a pointing error of 0.1° 
would produce a loss of signal ranging from 0.1% for N^^^ = 
20-39 to 6.0% for Npa^ >2000. As a consequence, the spec¬ 
tral index would increase by 0.01 and the flux normalization 
would decrease by 2%. 

3) Background evaluation. Our measurement is based on 
a very precise evaluation of the background. As explained 
in Section 3, the number of gamma rays is given by the dif¬ 
ference between the number of events detected in the event 
map (that contains the source events plus the cosmic ray back¬ 
ground) and the number of background events estimated with 
the time swapping method. Since the ratio between the num¬ 
ber of gamma rays and the number of background events is 


very small (ranging from ^3 x 10 for N^a^ = 20-39 up to 
^4 X 10“^ for Npa^ > 300), even a small systematic error 
in the background evaluation could produce a big error in the 
source flux. 

Possible sources of systematics are; a) the presence of cos¬ 
mic ray excess regi ons due to the medi um scale anisotropy, 
as those reported in (IBartoli et alJl2013'3l . close to the source, 
b) changes in atmospheric conditions able to modify the back¬ 
ground distribution in local coordinates in less than 2-3 hours, 
and c) the detector malfunctioning. 

Such effects, when present, could generate extended re¬ 
gions in the signal map with evident excesses or deficits, in 
some cases involving the whole map. Instead, an accurate 
evaluation of the background produces a map with all the bin 
contents consistent with zero except at the position of real 
sources. 

Concerning the medium scale anisotropy, we adopted a par¬ 
ticular procedure to correct the background systematics in the 
sky regions coincid ent or adjacent to cosmic ray excesses 
(IBartoli et alJl2013h . This correction is not necessary in the 
Crab Nebula region. 

Concerning points b and c, it has to be specified that the 
maps are built with data sets of 2-3 hours and individually 
checked. When a map shows significant anomalies, the cor¬ 
responding dataset is rejected, so only “good maps” are com¬ 
bined to build to the “total maps”. 

To test the background reliability of the total maps, we use 
the regions that are not involved in the Crab Nebula emission, 
i.e., the bins with an angular distance from the source larger 
than a minimum value, depending on the PSF. From these 
“out-source” regions we expect no significant excess, since 
they do not include any other known gamma ray source with 
a flux above the ARGO-YBJ sensitivity. For any of the nine 
maps we have evaluated the distribution of the excesses in the 
out-source region bins (before smoothing). We found that all 
the distributions are well described by Gauss functions with 
mean values consistent with zero and r.m.s. consistent with 
unit. 

Adding all the nine maps together, the total number of 
events detected from the out-source regions is 1.18 x 10®. 
This value differs from the corresponding estimated back¬ 
ground by -9.3 X 10^ events, corresponding to -0.3 standard 
deviations. Since there is no significant excess or deficit of 
events, we can calculate the upper limit of the systematic er¬ 
ror in the out-source region. We found that the relative error 
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in the background value is less than 3.7 x 10“^ at 90% confi¬ 
dence level. 

We can reasonably assume that a similar systematic error 
involves the region of the map containing the source signal, 
and that all the nine maps have comparable systematic errors. 
Based on these assumptions, we can evaluate the effects of 
such an error on the signal, which are obviously more rel¬ 
evant for the maps in which the signal-to-background ratio 
is smaller. We found that the error in the photon number 
is <13% for = 20-39, <1% for = 100-199, and 
<0.01% for > 1000. 

According to these values, the corresponding systematic er¬ 
ror in the spectrum flux normalization would be less than 2%, 
and the error in the spectral index would be less than 0.05. 

4) Event rate variations. Studying the rate of cosmic ray 
showers over five years, we observed variations on timescales 
from hours to months up to 10% with respect to the mean 
value. These variations are mostly due to; a) variation of at¬ 
mospheric pressure and temperature that modify the showers 
propagation in the atmosphere, b) variation of the detector ef¬ 
ficiency due to changes of the local temperature and pressure, 
c) aging of the detector. 

Gamma rays are assumed to be subject to similar variations. 
To study the stability of the Crab Nebula flux, we corrected the 
rate of the events observed from the source using the cosmic 
ray rate as a normalization factor (see Section 5.3). However, 
an absolute normalization cannot be performed. The Monte 
carlo simulations refer to a fixed atmospheric condition and 
a given detection efficiency, that cannot exactly reproduce 
the average effect over several years of different conditions. 
Considering the amount of the observed rate variation, a rea¬ 
sonable estimation indicates a possible systematic error in the 
flux smaller than 4%. 

Total systematic error. Adding all these contributions lin¬ 
early, we conservatively estimate the total systematic error to 
be less than 30% for the flux normalization and 0.06 for the 
spectral index. 


5. CRAB NEBULA LIGHT CURVE 

To Study the stability of the Crab Nebula emission, we con¬ 
sider the events with Npad > 40. Our total dataset consists of 
1851 days, with an average observation time of 6.0 hours per 
day. The average rate of events with Npaj > 40 is 137 ± 10 
day“'. 

Figure [T] shows the observed rate for events with N^a^ > 
40, 100 and 500, as a function of the Julian date, in bins of 
200 days. The median energies corresponding to these N^aj 
thresholds are 0.76, 1.8 and 5.1 TeV, respectively. The sig¬ 
nal appears stable during five years for any threshold within 
the statistical fluctuations. Assuming a constant rate, the ob¬ 
tained 15.7, 3.27 and 5.17 (with 9 d.o.f.) for N^^^ > 

40, 100 and 500. The corresponding p-values are 0.073, 0.95 
and 0.82, respectively. 

A six-year monitoring of the Crab Nebula was previously 
performed by the Tibet-Ill air shower array from 1999 to 
2005, at energies ~3 TeV, with a sensitivity 3-4 times lower 
than that of ARGO-Y BJ, reporting a yearly fl ux consistent 
with a steady emission (lAmenomori et al.ll2009l) . 

5.1. Search for Flares 

To make a “blind” search for short time rate variations, we 
consider all the time intervals of duration Af ranging from 1 to 
15 days, starting from every observation day. This time range 


5 200 

W 

I 150 
100 


50 

0 


40 

35 

30 

25 

20 

15 

10 

5 

0 

5 

4 


3 


2 
1 
0 

4400 4600 4800 5000 5200 5400 5600 5800 6000 6200 6400 

Time (MJD-50000) 

Figure 7. Rate of events detected from the Crab Nebula for different Np^^ 
thresholds as a function of time. The bin width is 200 days. The lines repre¬ 
sent the average values. 
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has been chosen on the basis of the duration of the flares ob¬ 
served in the GeV energy region. For each interval we com¬ 
pare the observed rate of Crab events with the average rate 
and evaluate the significance of the excess as: cr, = {Ri-R,n) 
/ S{Ri-R,n), where R, is the counting rate in the i-th interval, 
R,n is the average counting rate, and i5(R, -R,„) is the statistical 
error of the difference Ri-R,„. Note that the values of a, are 
not ind^endent, since the time intervals overlap. 

Fig. [8] shows the distributions of cr, for Af = 1 day and At 
= 2-15 days, for Npad > 40. The total number of intervals is 
1851 for Af = 1 and 25911 for Af = 2-15 days. The distri¬ 
butions can be fitted by a Gauss function with mean value m 
= -0.04 ± 0.03 and r.m.s.= 1.05 ± 0.02 for Af = 1 day, and 
m = -0.06 ± 0.01 and r.m.s.= 1.061 ± 0.005 for Af = 2-15 
days. The root mean square values indicate rate variations 
slightly larger than what expected by statistical fluctuations. 
However, no significant excess is observed for any of the con¬ 
sidered time intervals. 

Given the ARGO-YBJ sensitivity, a flare would produce a 
5 s.d. signal (pre-trial) if the flux exceeds the average value 
by a factor / ~ 10 / •y^Af(days). 

5.2. Correlation with Fermi-LATData 

To reduce the number of trials in the search for possible 
flares, we can limit our analysis to the days in which a flare 
was observed by satellite instruments at lower energies. We 
consider the Fermi-LAT daily light curve at energy E > 100 
MeV from 2008 August to 2013 February, obtained through 
the analysis of the scientific Fermi data publicly available at 
the Fermi Science Support Centefl 

The first panel of Figure |9] shows the daily light curve, rep¬ 
resenting the sum of the nebula and pulsar fluxes. The average 

* http://fermi.gsfc.nasa.gov/ssc/data/access/ 
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Significance (s.d.) 


Figure 8. Search for flares: distribution of the excesses from the Crab Nebula 
around the average value, in units of standard deviations, for different flare 
durations At. 



Time (MJD-SOOOO) 


Figure 9. Panel 1 : daily light curve of the Crab Nebula and pulsar by Fermi- 
LAT; Panel 2: light curve of the Crab Nebula (pulsar subtracted) in bins of 
200 days, by Fermi-LAT; Panel 3: light curve of the Crab Nebula by ARGO- 
YBJ for events with Npad > 40. The dashed line in panels 2 and 3 has been 
obtained excluding the days with flai'es. 

flux is (2.66 ± 0.01) X 10“^ photons cm“^ s“'. Also excluding 
the days with flares, the rate is variable, with modulations on 
timescales of weeks and months. 

First we consider the three largest Fermi flares, which oc¬ 
curred in 2009 Fe bruary, 2010 Septem ber and 2011 April 
dAbdo et alj|201 Ih iBuehiler et al.ll2012l) . To define the time 
boundaries and the duration of these flares we select the days 


in which the Fermi flux is higher than 4x10^ photons cm ^ 
s“'. The dates and the duration of the three flares are given in 
Table [3 The counting rates from the Crab Nebula measured 
by ARGO-YBJ with events with > 40 during the flares 
are compared with the average rate of 137 ± 10 day“'. In all 
cases, the rates are slightly higher than the average value, but 
consistent with it within statistical errors (see Table O. Sum¬ 
ming the three flares the average rate is 205 ±91 day“*. Table 
[3 also shows the results concerning the events with Npad > 
100 and 500. No significant excess is present in this case, 
either. 

Our preliminary analysis reported in (lAielli et alJ 1201 Ol) 
showed a 4 standard deviations excess observed in the time 
interval from 2010 September 17 to 22 from a direction con¬ 
sistent with the Crab Nebula. However, removing the contri¬ 
bution of the steady flux and taking into account the number of 
trials, the post-trial significance was about two standard devi¬ 
ations. A further excess with a similar post-trial stati stical sig¬ 
nifica nce was observed during the 2011 April flare (IVernettol 
l2013h . In the present work, based on a better shower recon¬ 
struction and the event selection described in Section 2.2, the 
significance of the Crab Nebula signal integrated over 5 years 
increases by about 15% with respect to the old analysis, but 
the signal observed during the Fermi flares decreases. The 
flux measured during the flares appears slightly higher than 
what would be expected from the steady emission, but consis¬ 
tent with it within one standard deviation. Both our previous 
analysis and the current one hint at a possible flux enhance¬ 
ment during the flares, but the reduced significance prevents 
us from drawing a definitive conclusion. 

To extend the search for flares to the whole observation 
time, and not limit the analysis to the largest flares, we se¬ 
lected the Fermi data according to the measured daily flux and 
checked the corresponding ARGO-YBJ event rate. Table |4] 
reports the ARGO-YBJ rates for different levels of the Fermi 
flux and different Npa^; thresholds. The rates are consistent 
with the average rate for any Fermi flux level. In particular, 
the ARGO-YBJ rate (for Npa^; > 40) measured in the 62 days 
in which the Fermi flux exceeds 4 x 10“® ph cm“^ s“', is 190 
± 55 day“* i.e. 1.4 ± 0.4 times higher than the average rate. 

Finally, to study a possible correlation on timescales of 
months or years, we compare the light curves of the two detec¬ 
tors over the common observing time (^ 4.5 years), dividing 
the data into bins of 200 days. The bin width is chosen in or¬ 
der to have a significant signal in the ARGO-YBJ data (about 
7 s.d.). 

Since the flux measured by Fermi is the sum of the nebula 
and pulsar contributions, and since the pulsar flux Fp aver¬ 
aged over the pulsation period is also stable during flares {Fp 
= (2.04 ± 0-00 X 1 0~^ photons cm“^ s“' for E >100 MeV 
(iBuehler et al.ll20f2h ). the flux of the pulsar has been sub¬ 
tracted. The obtained average nebula flux is (6.2 ±0.1) x 
10“^ photons cm“^ s“'. The nebula flux shows variations up 
to ^30% of the average flux, with = 126, for 8 d.o.f. (see 
the second panel of Fig.|3l- The large variations are not only 
due to flares. In the same figure the dashed curve shows the 
flux obtained excluding the 62 “flaring days”. In this case the 
average value is (5.6 ± 0.1) x 10“^ photons cm“^ s“', with 
= 80. 

The lower panel of Fig.|3shows the corresponding ARGO- 
YBJ data for Npa^; > 40. The average rate is 139.3 ± 10.6 
events day“' (x^ = 14.2 for 8 d.o.f., p-value p = 0.077). Even 
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Table 3 

ARGO-YBJ results during the three largest Crab Nebula flares observed by Fermi-LAT. 



At 

(days) 

Max.Fermi flux 
(ph cm”^ s“^) 

ARGO-YBJ 
Observation 
time (hr) 

ARGO rate 
>40 
(ph day“*) 

ARGO rate 
>100 
(ph day“^) 

ARGO rate 
Nparf >500 
(ph day“^) 

Flare D 

8 

6.3±0.8x10-6 

49.6 

142 ± 151 

21 ±28 

2.5 ± 4.6 

Flare 2*’ 

5 

6.4±0.8xl0-® 

31.5 

265 ± 190 

58 ± 36 

-3.8 ±5.7 

Flare 3“^ 

9 

19.8±0.8xl0-'S 

58.0 

228 ± 144 

51 ±27 

2.9 ± 4.4 

Sum of 3 flares 

22 


139 

205 ± 91 

41 ± 17 

1.3 ±2.8 

All ARGO data 




137 ± 10 

27 ±2 

2.8 ±0.3 


Note. — 

“ Start time MJD 54864 (2009 February 02) 
Start time MJD 55457 (2010 September 18) 
" Start time MJD 55662 (2011 April 11) 


Table 4 

ARGO-YBJ photon rate for different flux levels measured by Fermi. 


Fermi flux 
( 10^ phcm-2 s*') 

Number of days ARGO rate 
>40 
(ph day“‘) 

ARGO rate 
Np„rf>100 
(ph day“‘) 

ARGO rate 

>500 

(phday-*) 

< 2.0 

175 

197 ± 34 

26 ±6 

3.1 ± 1.0 

2.0-3.0 

915 

118 ± 15 

27 ±3 

2.8 ±0.4 

3.0-4.0 

435 

148 ± 21 

28 ±4 

3.0 ±4.3 

4.0-5.0 

46 

188 ±64 

33 ± 12 

1.9 ±2.0 

> 5.0 

16 

198 ± 107 

54 ±20 

0.2 ± 3.2 

All ARGO data 


137 ± 10 

27 ±2 

2.8 ±0.3 


if the ARGO-YBJ rate variations are consistent with statistical 
fluctuations, the Fermi and ARGO-YBJ data seems to follow 
a similar trend. The ARGO-YBJ rate appears higher in the 
“hot” Fermi periods. The dashed curve is obtained after the 
exclusion of the flaring days. 

Fig. [TOl shows the ARGO-YBJ percentage rate variation 
with respect to the mean value (AFarco) as a function of the 
corresponding variation of the Fermi rate (AFfermi), for the 
9 bins of the light curve. The Pearson correlation coefficient 
between the two data sets is r = 0.56 ± 0.22. The quoted error 
for r is the root mean square of the distribution of the correla¬ 
tion coefficients obtained by simulating the fluctuations of the 
counting rates of each bin, according to their statistical errors. 
Fitting the 9 points with the function AFargo = « AFfermi + b, 
the values of the best-fit parameters are a = 0.88 ± 0.37 and b 
= 0.018 ± 0.079, with = 8.3 for 7 d.o.f. Discarding the 62 
“flaring” days the correlation coefficient becomes r = 0.45 ± 
0.23 and the parameters of the linear fit are a = 0.96 ± 0.45 
and b = 0.018 ± 0.082, with = 10.4. 

The same analysis has been performed using a different bin 
width, ranging from 10 to 450 days. The corresponding cor¬ 
relation coefficient steadily increases from r = 0.10 ± 0.06 
(10 days) to r = 0.59 ± 0.23 (450 days). It has to be noted, 
however, that when using a small bin width, the ARGO-YBJ 
signal is not significant enough to search for a correlation un¬ 
less the flux variations are very large. In 10 days, for example, 
the average ARGO-YBJ signal is 137 ± 135 events day“'. 
The statistical fluctuations would hide a possible flux varia¬ 
tion, unless the flux becomes more than a factor of 4-5 higher 
than the average. 

The above results refer to events with > 40. The cor¬ 
relation coefficient is lower when selecting more energetic 


events: using a bin width of 200 days, for Fipad > 100, r = 
0.19 ± 0.31; and forN^^^ > 500, r = 0.46 ± 0.28. 

5.3. Stability of the ARGO-YBJ Data 

When studying the time evolution of a signal over several 
years, a discussion on the possible causes of detector insta¬ 
bilities is mandatory, to exclude systematic effects that could 
produce artificial rate variations. Since the measured number 
of events from the source Ns = Nf-N r is the difference be¬ 
tween the number of events Nf detected in the source map and 
the number of background events Nr estimated with the time 
swapping method, one must separately analyze the stability of 
the different contributions. 

1) A loss of signal events Ns could be produced by vari¬ 
ations of the pointing accuracy. Studying the Moon shadow 
month by mont h, we have verified that the pointing is stable 
within 0.1 deg (iBartoli et al.ll201 fbh . Given the moderate an¬ 
gular resolution for events with Npad > 40, such a value could 
produce signal fluctuations of less than 2%. 

2) A worsening of the detector angular resolution (due to an 
increase of the time resolution of RPCs occurring at particu¬ 
larly low temperatures) could produce a loss of signal events 
Ns- A broadening of the PSF would also cause a decrease of 
the Moon shadow signal. That, however, is found to be stable 
within statistical fluctuations. 

3) Atmospheric pressure and temperature variations can af¬ 
fect the RPC detection efficiency, which can also be altered 
by some RPC not working properly or by aging effects. 

4) Pressure and temperature produce changes in the shower 
rate of the order of a few percent due to the different condi¬ 
tions in which the showers propagate in the atmosphere. 

The two latter effects modify Ns, Nf and Nr by about the 
same factor (neglecting the different behavior of cosmic ray 
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Figure 10. Percentage variation of the Crab Nebula flux with respect to the 
average value: ARGO-YBJ vs Fermi-LAT data. The straight line is the best 
fit curve. 


and gamma ray showers, that in this contest can be considered 
a second-order effect). This allows the use of Nb to correct 
the Crab rate, multiplying the Crab rate observed in a given 
time interval by the correction factor fc = B^IB, where B^ is 
the average background rate and B is the background rate in 
that interval. The light curve in Figure [T] has been corrected 
according to this method, with fc ranging from 0.91 to 1.07. 

5) Further possible systematics could be an incorrect evalu¬ 
ation of the background In Section 4 we evaluated the ac¬ 
curacy of the background for the total source signal. To check 
the accuracy of the background along the years, we can use 
the same out-source regions previously defined. For events 
with Npflrf > 40, the out-source light curve in 200-day bins 
has a mean value of -7.9 ± 19.0 events day“* and a X' = 10.2 
for 9 d.o.f., corresponding to a p-value p = 0.67. According 
to these results the background is stable and should not intro¬ 
duce any systematic effect on the rate of the Crab signal. 

6. SUMMARY AND CONCLUSIONS 

The ARGO-YBJ events recorded over five years have been 
analyzed to evaluate the Crab Nebula spectrum and study 
the temporal behavior of the gamma ray emission. Using 
the events with '^pad > 20, the statistical significance of the 
gamma ray signal reaches more than 21 standard deviations, 
and the observed photon rate is 189 ± 16 day”'. The event 
angular distributions around the source are well described by 
the PSFs obtained by simulations. 

The source spectrum extends over nearly 2 decades in en¬ 
ergy and five decades in flux. The spectral shape is consistent 
with a power law behavior in the range 0.3-20 TeV with a 
spectral index -2.63 ± 0.05. An exponential cutoff would be 
consistent with our data in case of a cutoff energy higher than 
12 TeV at 90% confidence level. 

The study of the Crab Nebula light curve has been carried 
out to check the stability of the flux over years and to search 
for possible flares on the timescale of days. All the known 
sources of rate instabilities have been examined and the ef¬ 
fects corrected. 


Concerning flares, a blind search for flux increases of dura¬ 
tion between 1 and 15 days shows no significant excess. The 
average rate of events with N^a^ > 40 measured by ARGO- 
YBJ during the three most powerful flares detected by Fermi- 
LAT (in 2009 February, 2010 September and 2011 April) is 
205 ±91 day”*, which is consistent with the average value of 
137 ± 10 day”'. 

The five year ARGO-YBJ light curve with a binning of 200 
days is consistent with a constant flux with a probability of 
0.07. A correlation analysis with the corresponding Fermi- 
LAT data gives a Pearson correlation coefficient r = 0.56 ± 
0.22. The small statistical significance of these results does 
not allow the claim for a flux variability correlated with the 
observations at lower energies. If such a correlation was due 
to a real astrophysical phenomenon, the found regression co¬ 
efficient a = 0.88 ± 0.37 would imply a similar percentage 
variation in Fermi and ARGO-YBJ rates, suggesting a similar 
behavior of the gamma ray emission at energies ^100 MeV 
and-1 TeV. 

So far, no variation of the Crab Nebula flux at TeV ener¬ 
gies has been reported by any detector. Assuming the flares 
observed by AGILE and Fermi due to synchrotron radiation 
from a population of electrons accelerated up to 10*^ eV, 
the Inverse Compton emission associated with this population 
would occur in the Klein-Nishina regime and would produce 
gamma rays of energy approximately equal to that of the elec¬ 
trons. Such a flux would not be detectable by any of the exist¬ 
ing gamma ray experiments. With these assumptions, a TeV 
excess could hardly be intepreted as IC emission associated 
with the synchrotron radiation observed at lower energies, and 
would require a completely new interpretation. 
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